Revised, September 20, 2000 

Preprint typeset using LATj^X style cmulateapj v. 20/04/00 



o- 
o: 

O. 
(N ■ 

<D ■ 
C/V 



> 

o: 

^ ■ 

OV 

o: 
o. 

Of 

i ' 

o- 

a: 



13 



WATER AND DUST EMISSION FROM W HYDRAE 
Victor Zubko 1 and Moshe Elitzur 

Department of Physics and Astronomy, University of Kentucky, Lexington, KY 40506-0055 

Revised, September 20, 2000 

ABSTRACT 

We construct a self-consistent model for the wind around W Hya by solving the coupled equations 
describing the hydrodynamics and dust radiative transfer problems. The model matches simultaneously 
the observed continuum radiation and wind velocity profile. The water line emission is calculated next 
using the water abundance as the only free parameter, fitted from the ISO observations of Neufeld et al. 
(1996) and Barlow et al. (1996). The gas temperature is determined from a thermal balance calculation 
that includes water as one of its main components. Our model successfully fits all the observed water 
lines, resolving a major discrepancy between the modeling results of the two observing teams. The mass 
loss rate is 2.3xl0 -6 M Q yr _1 , the water abundance is l.OxlO -4 and the ortho:para ratio is 1:1.3. 

Subject headings: stars: AGB and post-AGB — stars: individual: W Hya — infrared: stars — 
circumstellar matter — dust, extinction — molecular processes 



1. INTRODUCTION 

Water is a dominant coolant of outflows around cool 
oxygen-rich stars (Goldreich & Scoville 1976, GS hereafter; 
Chen & Neufeld 1995, CN hereafter; Truong-Bach et al 
1999, TB hereafter). However, until recently observations 
of water cooling lines were impossible because of the at- 
mospheric opacity at these wavelengths. The situation has 
changed with the successful launch of the Infrared Space 
Observatory (ISO). One of the first objects observed with 
ISO was W Hydrae, an M7.5 semi-regular red giant, and 
it showed the expected thermal water emission in obser- 
vations with both the SWS (Neufeld et al 1996; hereafter 
N96) and LWS instruments (Barlow et al 1996; hereafter 
B96). Both teams also fitted their observations based on 
the GS approach, resulting in strikingly different estimates 
for the mass loss rate: 6x10 M Q yr _1 (B96) and (0.5- 
3)xl0 -5 M Q yr _1 (N96). Here we aim to resolve this 
discrepancy. In contrast with the original studies we con- 
struct a self-consistent model of the radiation field, dust, 
and gas in the shell, taking account of all the infrared con- 
tinuum observations as an additional constraint. 

2. MODELING 

The driving force of the wind is radiation pressure on 
the dust, the gas particles are dragged along by collisions 
with the dust grains. The internal properties of the gas, 
such as temperature, do not play any role in the dynamics, 
the wind structure can be obtained by solving the coupled 
equations for hydrodynamics and dust radiative transfer. 
We now describe our calculation for W Hya. With the de- 
rived model we proceed to solve the H2O level population 
problem. 



2.1. 



Dynamics and IR Emission 

A complete calculation of the wind structure requires 
a solution of the coupled hydrodynamics and dust radia- 
tive transfer problems. Traditionally these calculations 
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involved a large number of input parameters. However, 
Ivezic & Elitzur (1995) noted that the dusty wind prob- 
lem possesses general scaling properties such that, for a 
given type of grains, both the dynamics and radiative 
transfer depend primarily on a single parameter - the 
overall optical depth. Subsequent analysis by Ivezic & 
Elitzur (1997) established rigorously that the dust radia- 
tive transfer problem possesses scaling properties under 
the most general circumstances. Scaling was incorporated 
in the code DUSTY 2 (Ivezic et al. 1999), which solves fully 
the dusty wind problem. The solution provides the radial 
variation of the velocity and radiation fields in terms of 
the scaled distance y — r/R- ln , where R ln is the shell in- 
ner boundary. That boundary is defined by the condition 
T c = Tdust(y = l)j where T c is the dust condensation tem- 
perature. The actual value of R m never enters. 

We use for modeling "astronomical silicate" dust grains 
with optical constants from Laor & Draine (1993) and the 
power law size distribution of Mathis et al. (1977). We 
assume prompt dust formation at y = 1 and no further 
grain growth or destruction. The only input parameter in 
addition to the dust properties is the stellar temperature 
T* = 2500K (Haniffet al. 1995). From a series of DUSTY 
models in which we varied the visual optical depth ry, the 
temperature T c and the shell outer radius Y — R ut/Rm 
we chose the one that best fits all the observations all the 
way from 1 /iin to 1.2 mm; because of its irregular vari- 
ability, optical data is not included. We find that models 
with r v = 0.7-1.0, T c = 900-1100 K and Y > 1000 are 
almost equally successful in reproducing the observational 
spectrum. The most significant parameter by far is tv; T c 
has only a small effect on the fitting and the role of Y is 
marginal. The model with r v = 0.83, T c =1000 K and Y = 
12,000 minimizes the fitting errors and is presented in Fig. 
1. Matching the model flux to observations in scale as well 
as spectral shape determines the shell angular dimensions 
(Ivezic & Elitzur 1997). The shell inner diameter is 0'.'182, 
consistent with the measured stellar diameter 0'.'046 (Han- 
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iff ct al. 1995), and its outer diameter 36', in agreement 
with the IRAS observations of Hawkins (1990). 
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Fig. 1. — The spectral energy distribution of W Hya. Data in- 
dicated by filled squares (NIR) are from Wilson et al. (1972); open 
squares (IRAS) from Hawkins (1990); filled circles (JCMT) from 
van der Veen ct al. (1995) and Walmsley et al. (1991). The thin 
solid line shows the ISO data of Barlow ct al. (1996). The bold 
solid line is our model result. 

2.2. Water Lines 

The statistical rate equations for the water level popu- 
lations require the following input: the radiation field, to 
determine radiative excitations; the gas density and tem- 
perature, for the collision rates; and the water column den- 
sity, required for the line radiative transfer (Elitzur 1992). 
DUSTY's output provides the radiation field as well as 
the dimcnsionless velocity and density profiles in the shell. 
Fixing the scale of densities requires two additional input 
properties. We take the distance d = 115 pc from the Hip- 
parcos catalog. This sets i?; n = 1.6xl0 14 cm and R out ~ 
1 pc, and the luminosity is 11,050 L Q in agreement with 
Haniff et al. (1995). Next, the wind terminal velocity v c = 
8 km s _1 (Young 1995) fixes the velocity scale and deter- 
mines the gas-to-dust mass ratio r g d = 850 and the mass- 
loss rate M = 2.3 xl0~ 6 M Q yr _1 . The complete velocity 
profile is shown in figure 2 together with the data for OH 
and H2O masers as well as CO thermal emission. All are 
properly explained by the model results. The SiO data is 
displaced from the wind velocity profile, as expected for 
this maser's location inside the dust formation zone (cf 
Elitzur 1992). 

The only required input quantities that remain unknown 
are the gas temperature and water abundance, and we 
fit those simultaneously from the water line observations. 
We calculate the temperature from the balance of cooling 
and heating due to adiabatic expansion, grain-gas colli- 
sions and H2O rovibrational transitions (we estimate the 
H2 vibrational contribution and find it negligible). We 
solve for the populations of the lowest 45 rotational lev- 
els of the ground vibrational state of ortho- and para- 



water; this accounts for all levels with energy < 2000 K 
above ground. The molecular data are from the HITRAN 
database (Rothman et al. 1998), the collision rate coeffi- 
cients from Green et al. (1993). The level populations and 
line emissivities are calculated with the escape probabil- 
ity method as functions of distance r, and the line fluxes 
by integrating the line emissivities over the shell volume. 
The free parameters are the abundances of the two wa- 
ter species, which must be considered independently since 
there is no radiative coupling between them. Detailed 
modeling shows that these abundances are constant in the 
region where the water emission originates (CN, TB). Be- 
cause of the central role of H2O in the energy balance, 
the calculations of the temperature and water line emis- 
sion are coupled, and repeated until the best fit achieved 
for the line observations. The best fit is found for water 
abundance n(H20)/n(H2) = 1.0xl0~ 4 and the ratio or- 
tho:para = 1:1.3. Figure 3 shows the temperature profile. 
The model parameters are summarized in Table 1, results 
and comparison with observations in Table 2. 
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Fig. 2. — The model prediction for the gas velocity profile. Data 
points (from Szymczak et al, 1998) show the H2O and OH mascr 
and CO thermal emission from the wind. The SiO maser emission 
originates in the extended atmosphere. 



3. DISCUSSION 

Our model fits all the water lines within the observa- 
tional errors, which generally exceed 50%. The quality of 
the fit for the LWS data is comparable to that of the B96 
model. For the SWS data N96 present a range of models 
for each line, and our model fits that data set better than 
any single one of them. Therefore, our model resolves the 
conflict among the previous water line calculations, fitting 
all the data with a single value for the mass loss rate. It 
is important to note that M is determined by the infrared 
data and v c , and remains unchanged during the model- 
ing of the water lines. We find that the acceptable range 
of M is ~ 2-3xl0~ 6 M Q yr _1 , with a nominal value of 
2.3xl0~ 6 M Q yr" 1 . Except for the high-end of the N96 
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range, most estimates of M are in agreement with ours 
within the errors (cf N96 and references therein). 
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Fig. 3. — The gas temperature profile from our best-fit model. 
Shown are also power-law fits for the Goldrcich & Scovillc (GS) and 
Chen & Neufeld (CN) temperature profiles. 

The source of the large discrepancy with N96 is not 
clear. Neufeld et al suggest that the temperature profile 
could be the reason, but this does not seem to be the case. 
Our temperature profile is not that different from the CN 
profile, which was used in the N96 study, especially in the 
relevant range T > 400 K (figure 3). We compared the 
various contributions to heating and cooling with those 
listed by CN and Neufeld & Kaufman (1993), and find 
good agreement; the difference in resulting profiles can be 
attributed to the different parameters used in the two cal- 
culations (CN used M = 3xl0~ 5 M yr _1 ). The two 
profiles differ much more with the GS profile employed in 
the B96 model, which produced M similar to ours. We 
suspect that a more important source of difference could 
be the radiation field, since radiative excitations play an 
important role in the water population distribution. A 
comparison is impossible since N96 do not give details of 
the radiation field they employed. A proper radiation field 
appears crucial for the water line calculations. 

The uncertainty in fitted parameters depends on the 
particular observations that constrain them. We estimate 
that the acceptable range for r gc j is ~ 850±100. We exper- 
imented with both power-law and single-size (a = O.f /jm) 
grain distributions and found the differences negligible. 
The large uncertainties in the measured line fluxes trans- 
late into a large uncertainty in the abundances of the two 
species of water. Between the two, the para-ri^O abun- 
dance is subject to the larger uncertainty and can vary by 
as much as factor 3, so that the ortho:para ratio can be 
anywhere from f to f /3. By comparison, B96 obtained 1 
for this ratio. All of these results differ greatly from the 
thermodynamic limit of 3. Depending on the ortho:para 
ratio, we find acceptable models for n(H20)/n(H2) in the 
range ~ 1-4 xl0~ 4 . 

In the modeling efforts of both B96 and N96, M was 
one of numerous free parameters fitted from the water line 



observations. The scaling approach taken here reduces 
the number of free parameters to the essential minimum 
and determines M prior to the water line fitting. Table 1 
breaks the model parameters into three categories. Quan- 
tities in the first division are specified as input, in addition 
to the IR and water observations. The second group in- 
cludes the free parameters of our two fitting procedures 
while the third group lists results derived from our fits. 
Given grain properties, the IR observations are fitted with 
the three free parameters tv, T c and Y and the single in- 
dependent input T*. The resulting model determines also 
the dimcnsionless velocity profile. Adding as independent 
input the source distance and the wind final velocity, the 
model results determine the full velocity profile, and Ty de- 
termines also M and the gas-to-dust ratio. Fitting the wa- 
ter lines involves no other input and only two free parame- 
ters — the ortho- and para-water abundances. Thanks to 
the central role of water lines in the gas temperature calcu- 
lation, the gas temperature is determined self-consistently 
as part of this second fitting procedure. Since the only 
free parameters in fitting the observed water fluxes are 
the abundances of the two species, confidence in the de- 
rived values is greatly enhanced. In principle, we could 
have used M as a free parameter in the water calcula- 
tions, as in the previous studies. In that case, consistency 
between the results of the two fitting procedures would be 
used as an additional constraint, automatically met by our 
calculation. 

While radiation pressure on the dust grains is generally 
accepted as the driving mechanism behind the wind ex- 
pansion, this mechanism has not been fully tested. The 
solutions of radiative transfer and the hydrodynamics 
problems must result in the same parameters, but this 
fundamental test has not been performed thus far; self- 
consistent modeling of both the IR emission and the wind 
structure in the same source has not yet been attempted. 
The most detailed previous calculation we are aware off 
is the TB modeling of R Cas. However, in that work the 
radiation field was not calculated self-consistently, instead 
it was derived from a dust temperature profile that was 
assumed beforehand as an input property. In contrast, 
DUSTY determines this temperature from a proper cal- 
culation of radiative equilibrium coupled to the radiative 
transfer including dust scattering, absorption and emis- 
sion. The spectral energy distribution is fitted with just 
three free parameters, only one of which (tv) is signifi- 
cant. Once these parameters are set, the outflow terminal 
velocity determines the entire velocity profile without any 
more freedom in the model. It is highly significant that 
a single self-consistent model with the minimal necessary 
number of parameters provides agreement with both the 
spectral energy distribution and the molecular velocity ob- 
servations. Apart from resolving the M discrepancy in W 
Hya, the success of our model provides strong support for 
the basic paradigm of winds in late-type stars. 

We thank Michael Barlow for the ISO LWS spectrum 
of W Hya, David Neufeld for help with the water molecu- 
lar data, Zeljko Ivezic for his comments on the manuscript 
and the referee for useful suggestions and for pointing out 
the Hipparcos distance measurement. The partial support 
of NASA and NSF is gratefully acknowledged. 
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Table 1 
Model Parameters 



Input parameters: 






Stellar temperature 3. 


2500 K 




Distance 15 


115 pc 




Final velocity 


8 km s _1 




Fitting parameters: 






Optical depth at 0.55 /zm 1 


0.83 




Dust temperature at R[ n l 


1000 K 




Shell thickness {Rout/Rin) 1 


12,000 




H2O abundance at i?i n w 


1.0xl0~ 4 




ortho-H 2 0:para-H 2 O w 


1:1.3 




Derived parameters: 






Stellar radius 


3.9xl0 13 cm 




Stellar luminosity 


11,050 L Q 




Shell inner radius (i?i n ) 


1.6xl0 14 cm 




Mass loss rate 


2.3 xl0~ 6 M Q yr" 


-1 


Gas-to-dust mass ratio 


850 




Shell outer radius (i? ou t) 


1.2 pc 




Dust temperature at i? ou t 


17 K 





a Haniff et al. (1995) 

b Hipparcos Catalog 

c Young (1995) 

1 fitted from IR observations 

w fitted from H 2 lines 
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Table 2 
Observed and model fluxes of H2O lines 



A (/im) 


Transition 


771 a 

r mod 


771 a 
-fobs 


F mnA 
F obs 


180.486 


0: 


2 2 i- 


->2 12 


3.43 


2.90 


1.18 


179.527 


0: 


2l2" 


-loi 


12.48 


8.66 


1.44 


174.624 


0: 


3o3~ 


->2 12 


8.31 


9.21 


0.90 


156.265 


0: 


523~ 


^4 32 


2.42 


5.41 


0.45 


144.517 


p: 


4i 3 - 


^3 22 


5.57 


6.38 


0.87 


132.408 


O: 


4 23 - 


->4i4 


2.97 


5.50 


0.54 


125.356 


p: 


4 04 - 


^3i 3 


12.39 


16.60 


0.75 


113.538 


O: 


4l4- 


^3o3 


11.27 


17.40 


0.65 


108.073 


0: 


2 2 i- 


-►llo 


15.47 


13.00 


1.19 


89.989 


p: 


322" 


-2n 


18.11 


27.30 


0.66 


83.283 


p: 


606- 


^5i 5 


14.42 


22.60 


0.64 


78.742 


0: 


4 23 - 


^3i2 


14.08 


28.00 


0.50 


67.089 


p: 


3 3 i- 


^2 2 o 


17.80 


59.80 


0.30 


66.438 


0: 


330- 


+2 2 1 


18.45 


22.90 


0.81 


63.458 


p: 


808- 


^n 


12.99 


30.40 


0.43 


58.699 


O: 


4 3 2" 


+3 2 1 


15.35 


19.50 


0.79 


57.637 


p: 


4 22 - 


+3 13 


30.13 


36.20 


0.83 


40.691 


O: 


4 32 - 


^3o3 


36.07 


23.00 


1.57 


31.772 


0: 


441- 


^3i2 


39.95 


63.00 


0.63 


29.837 


0: 


725- 


->6i 6 


30.72 


32.00 


0.96 


37.984 


0: 


441- 


■+4i4 


13.19 


28.00 


0.47 b 



Note. — Ortho transitions are marked by "o" , 
para by "p" . The 57-180 /im lines are from Bar- 
low et al. (1996), the 29-41 ^m lines from Neufcld 
et al. (1996). Q Fluxes are in 10~ 20 W cm~ 2 . 
''This line is severely contaminated by blending 
with two others. 



